Five different tree species were used to select parameter sets to represent different types found around Santa Barbara.
Steps taken to select the parameters were as follows:
## summarized annually
dfyr %>%
ggplot() + geom_line(aes(x=as.factor(year), y=lai, col=code, group=code))
dfyr %>% dplyr::filter(year>2009 & year < 2020) %>%
ggplot() + geom_line(aes(x=as.factor(year), y=lai, col=code, group=code))
dfyr %>%
ggplot() + geom_line(aes(x=as.factor(year), y=plantc, col=code, group=code))
#dfday %>% ggplot() + geom_line(aes(x=yd, y=lai, col=code))
dfyr2 %>% ungroup() %>%
dplyr::filter(year ==2011 | year==2014 | year==2017) %>%
ggplot() + geom_boxplot(aes(x=code, y=lai, fill=as.factor(year))) + ggtitle("annual avg LAI")
dfyr2 %>% ungroup() %>%
dplyr::filter(year ==2011 | year==2014 | year==2017) %>%
ggplot() + geom_boxplot(aes(x=code, y=plantc, fill=as.factor(year))) + ggtitle("annual avg plantc")
# load in mikes data
rs_dat <- read.csv("~/Google Drive/My Drive/plant_params/data/alonzo_data_df.csv")
rs_dat <- dplyr::filter(rs_dat, fractional_cov > 0.7)
trees_rs <- rs_dat %>% dplyr::filter(
class==81 | #QUAG
class==74 | #PLRA
class==73 | #PIUN
class==66 | #PICA
class==34 ) #EUGL
ggplot(trees_rs) + geom_boxplot(aes(x=as.factor(class), y=lai, fill=as.factor(class))) + ggtitle("LAI estimates from alonzo 2016")
## Warning: Removed 2143 rows containing non-finite values (stat_boxplot).
trees_rs$code <- NA
trees_rs$code[trees_rs$class==81] <- 'quag'
trees_rs$code[trees_rs$class==74] <- 'plra'
trees_rs$code[trees_rs$class==73] <- 'piun'
trees_rs$code[trees_rs$class==66] <- 'pica'
trees_rs$code[trees_rs$class==34] <- 'eugl'
ggplot() +
geom_boxplot(data=dfyr2, aes(x="model", y=lai, fill='rhessys'),position=position_dodge2(padding=0.5)) + facet_wrap('code') +
geom_boxplot(data=trees_rs, aes(x="data", y=lai, fill='alonzo'), position=position_dodge2(padding = 0.5)) + facet_wrap('code') +
ggtitle("compare LAI")
## Warning: Removed 2143 rows containing non-finite values (stat_boxplot).
ggplot() +
geom_boxplot(data=dfyr2, aes(x="model", y=plantc, fill='rhessys')) + facet_wrap('code') +
geom_boxplot(data=trees_rs[trees_rs$carbon<100,], aes(x="data", y=carbon, fill='alonzo')) + facet_wrap('code') +
ggtitle("compare plant carbon (kgC")
## Warning: Removed 2143 rows containing non-finite values (stat_boxplot).
# ggplot() +
# geom_boxplot(data=dfyr2, aes(x="model", y=lai, fill='rhessys')) + facet_wrap('code') +
# geom_boxplot(data=trees_rs, aes(x="data", y=lai, fill='alonzo')) + facet_wrap('code') +
# ggtitle("compare height (m)")
# load in David's data
vis <- read.csv("/Volumes/GoogleDrive/My Drive/plant_params/DM_data/tree_and_turfgrass_monroe2_selected_weightedmeanVIs.csv")
trees_vi <- vis %>% dplyr::filter(sp_class==81|
sp_class==74|
sp_class==73|
sp_class==66|
sp_class==34)
vi_dens <- ggplot(trees_vi) + geom_density(aes(x=NDVI, col=as.factor(year))) + facet_wrap("sp_code")
rh_dens <- dfyr2 %>% ungroup() %>%
dplyr::filter(year ==2011 | year==2014 | year==2017) %>%
ggplot() + geom_density(aes(x=lai, col=as.factor(year))) + facet_wrap("code")
grid.arrange(vi_dens, rh_dens)
Llorett et al. 2011 defines resilience as the capacity to reach pre-episode growth levels and relative resilience as the resilience weighted by the damage incurred during the episode.
By looking at time series as each year being relative the the pre-drought conditions, we can see how many years it takes to recover (or not at all). can also figure out which year to consider the drought occurring - where do we assign drought vs. post drought?
\(resistance = Drought / PreDrought\)
\(recovery = PostDrought / Drought\)
\(resilience = PostDrought / PreDrought\)
\(RelativeResilience = (PostDrought-Drought)/PreDrought\)
In the model, tree LAI responds the worst to the drought in 2015, after the months of experiencing the extreme heat in 2014 and 2015.
p_dec <- read.csv("decid_params.csv")
p_ev <- read.csv("evg_broadleaf_params.csv")
p_con <- read.csv("pica_conifer_params.csv")
p_dec2 <- cbind(p_dec[,15:54], p_dec$code)
p_ev2 <- cbind(p_ev[,15:54], p_ev$code)
p_con2 <- cbind(p_con[,15:54], p_con$code)
p_names <- c(
"pore_size_index",
"psi_air",
"run",
"epc.leaf_cn",
"epc.branch_turnover",
"epc.gl_smax",
"epc.flnr_age_mult",
"epc.litter_moist_coef",
"epc.psi_close",
"mrc.q10",
"epc.vpd_close",
"epc.leaf_turnover",
"epc.froot_turnover",
"epc.storage_transfer_prop",
"epc.height_to_stem_coef",
"epc.waring_pa",
"epc.root_distrib_parm",
"epc.proj_sla",
"epc.alloc_stemc_leafc",
"epc.ext_coef",
"specific_rain_capacity",
"epc.flnr_sunlit",
"epc.flnr_shade",
"epc.netpabs_sunlit",
"epc.netpabs_shade",
"epc.netpabs_age_mult",
"epc.cpool_mort_fract",
"epc.min_percent_leafg",
"epc.livewood_turnover",
"epc.waring_pb",
"epc.max_storage_percent",
"epc.min_leaf_carbon",
"epc.resprout_leaf_carbon",
"epc.alloc_frootc_leafc",
"epc.frootc_crootc",
"epc.alloc_prop_day_growth",
"epc.day_leafon",
"epc.day_leafoff",
"epc.ndays_expand",
"epc.ndays_litfall", "code")
colnames(p_dec2) <- p_names
colnames(p_ev2) <- p_names
colnames(p_con2) <- p_names
all_p <- rbind(p_dec2, p_ev2, p_con2)
###### now can use to plot
ggplot(all_p) + geom_boxplot(aes(x=as.factor(code), y=epc.proj_sla, fill=as.factor(code))) + ggtitle("proj_sla")
ggplot(all_p) + geom_boxplot(aes(x=as.factor(code), y=epc.leaf_cn, fill=as.factor(code))) + ggtitle("leaf_cn")
ggplot(all_p) + geom_boxplot(aes(x=as.factor(code), y=epc.waring_pa, fill=as.factor(code))) + ggtitle("waring_pa")
ggplot(all_p) + geom_boxplot(aes(x=as.factor(code), y=epc.branch_turnover, fill=as.factor(code))) + ggtitle("branch_turnover")
ggplot(all_p) + geom_boxplot(aes(x=as.factor(code), y=epc.leaf_turnover, fill=as.factor(code))) + ggtitle("leaf_turnover")
comb <- left_join(dfyr2, all_p, by=c('code','run')) %>%
dplyr::filter(year==2011 | year==2014 | year==2017) %>%
group_by(year, run, code) %>% summarize_all(funs(mean), na.rm=T)
key parameters for lai and plantc:
ggplot(comb) + geom_jitter(aes(x=epc.branch_turnover, y=plantc, col=as.factor(year))) + facet_wrap('code') +
ggtitle("branch_turnover")
ggplot(comb) + geom_jitter(aes(x=epc.leaf_turnover, y=lai, col=as.factor(year))) + facet_wrap('code') +
ggtitle("leaf_turnover") + geom_smooth(aes(x=epc.leaf_turnover, y=lai, col=as.factor(year)))
## `geom_smooth()` using method = 'gam' and formula 'y ~ s(x, bs = "cs")'
ggplot(comb) + geom_jitter(aes(x=epc.waring_pa, y=lai, col=as.factor(year))) + facet_wrap('code') + ggtitle("waring pa")
ggplot(comb) + geom_jitter(aes(x=epc.proj_sla, y=lai, col=as.factor(year))) + facet_wrap('code') + ggtitle("proj_sla")
ggplot(comb) + geom_jitter(aes(x=epc.psi_close, y=lai, col=as.factor(year))) + facet_wrap('code') +
ggtitle("psi_close")